This report examines the impact of scabies outbreaks on antiparasitic treatment trends and the preferred treatment options across Scotland during the second half of 2023 (June to December). The analysis focuses on prescription data to identify patterns and shifts in medication choices among health boards during this period.
To begin, the necessary R packages are loaded. These packages facilitate data cleaning, manipulation, visualization, and analysis. Following this, the monthly datasets for the study period are loaded. A vector of URLs corresponding to the CSV files for each month is created, ensuring that the process remains reproducible. The map() function is then used to iterate over these URLs, downloading the datasets, cleaning the column names for consistency, and storing the data in a list. Finally, the list elements are named according to the corresponding months, enabling straightforward access and identification of each dataset. them in a list. Finally, we assigned the names of the month that correspond to each dataset.
pacman::p_load(tidyverse, janitor, gt, gtExtras, purrr, ggrepel, sf, patchwork, here, leaflet, htmlwidgets, htmltools)
options(timeout = 1000)
urls <- c("https://www.opendata.nhs.scot/dataset/84393984-14e9-4b0d-a797-b288db64d088/resource/66e9611a-fbea-45b2-9edd-351c388fd06d/download/pitc202306.csv","https://www.opendata.nhs.scot/dataset/84393984-14e9-4b0d-a797-b288db64d088/resource/7bb45ee6-6f1c-45f4-933b-958bdbe0ca4f/download/pitc202307.csv","https://www.opendata.nhs.scot/dataset/84393984-14e9-4b0d-a797-b288db64d088/resource/72ca4ad2-0228-4672-9eb0-cc911e4a8ca7/download/pitc202308.csv","https://www.opendata.nhs.scot/dataset/84393984-14e9-4b0d-a797-b288db64d088/resource/2d3d240f-1467-4a91-9f0e-769745650cb9/download/pitc202309.csv","https://www.opendata.nhs.scot/dataset/84393984-14e9-4b0d-a797-b288db64d088/resource/94134438-db1d-4b8d-8f70-5e9b0e47bd03/download/pitc202310.csv","https://www.opendata.nhs.scot/dataset/84393984-14e9-4b0d-a797-b288db64d088/resource/21d11fed-a494-4e30-9bb6-46143c3b9530/download/pitc202311.csv","https://www.opendata.nhs.scot/dataset/84393984-14e9-4b0d-a797-b288db64d088/resource/00cdf8d7-f784-4e87-894c-34f2540ea6ab/download/pitc202312.csv")
data_list <- map(urls, ~ read.csv(.x) %>% clean_names())
dataset_year23 <- setNames(data_list, c("June", "July", "August", "September", "October", "November", "December"))
Additional datasets required for the analysis are loaded to support the study. These datasets include information on Health Boards, population data, GP practices, and spatial data. A list of Health Boards with their corresponding names is loaded, retaining only the Health Board code and name. Population estimates for each Health Board in 2023 are then loaded, and since the dataset is aggregated, only rows labeled “All” are selected. A dataset containing GP practice codes and their associated data zones is also imported. Spatial data defining the boundaries of Scottish data zones is loaded to provide geographic context. Finally, data zones are matched to Health Boards using a supplementary dataset.
nhs_board <- read.csv(here("data", "Health_Boards_(Dec_2020)_Names_and_Codes_in_Scotland.csv")) %>% clean_names()%>%select(hb = hb20cd, hb_name = hb20nm)
population_data2023 <- read.csv(here("data","hb2019_pop_est_14102024 (1).csv")) %>%clean_names() %>% filter(year == "2023" & sex == "All") %>% select(hb, all_ages)
gp_scot <- read.csv(here("data", "/practice_contactdetails_oct2023-open-data (1).csv"))%>% clean_names() %>% select(practice_code, data_zone, hb)
data_zones <- st_read(here("data", "SG_DataZone_Bdry_2011.shp"), quiet = TRUE)
DataZone_with_hb <-read_csv(here("data","DataZone2011lookup_2024-07-04 (1).csv")) %>%clean_names()%>% select(dz2011_code, hb_code)
A factor was created to include various treatments for scabies, categorized into topical (permethrin cream, malathion lotion, benzyl benzoate lotion, crotamiton cream, and sulfur ointment) and oral (ivermectin) options (British Association of Dermatologists, 2023). This factor was used to filter and analyze scabies treatment trends.
To analyze scabies treatment trends during the second half of 2023 (June to December), the monthly datasets were combined into a single comprehensive dataset. This allowed for a unified analysis of scabies-related prescriptions across health boards during the specified period.
The monthly datasets were filtered to include only prescriptions related to the defined scabies treatments. This step grouped the data by health board and GP practice, allowing the total quantity of medications prescribed each month to be calculated.Next, the processed dataset was joined with the Health Boards dataset to incorporate the corresponding board names. Subsequently, the population data for each health board was integrated into the dataset. This allowed for the calculation of the prescription rate per 1,000 people, enabling normalization of the data and ensuring comparability across different population sizes.Finally, the data was visualized through a line graph that displayed trends in scabies treatment prescriptions by health board over the specified months. Each health board was represented by a distinct line, with labels highlighting the values for the final month. This comprehensive approach ensures that prescription trends can be effectively compared across health boards while accounting for population differences.
scabies <- c("PERMETHRIN 5|PERMETHRIN_CR|MALATHI|BENZYL BENZOA|IVERMECTIN 3")
total_data23_sc <- dataset_year23 %>%
imap_dfr(~ .x %>%filter(str_detect(bnf_item_description, scabies) & !is.na(bnf_item_description)) %>% group_by(hbt, gp_practice) %>%summarise(paid_quantity = sum(paid_quantity, na.rm = TRUE), .groups = "drop") %>%mutate(month = .y))
with_names_sc <- total_data23_sc %>% full_join(nhs_board, join_by(hbt == hb)) %>% mutate(month =factor(month, levels = c("June", "July", "August", "September", "October", "November", "December"),ordered = TRUE)) %>%arrange(hb_name, month)
no_gp <- with_names_sc %>% filter(!is.na(month)) %>% group_by(month, hb_name, hbt) %>% summarise(paid_quantity = sum(paid_quantity))
rate_sc <- no_gp %>% full_join(population_data2023, join_by(hbt == hb)) %>% mutate(rate_per1k = paid_quantity/all_ages * 1000)%>%filter(!is.na(month), !is.na(hb_name))
sc_squares <- rate_sc%>% ggplot(aes(x = month, y = rate_per1k, colour = hb_name))+ geom_line(aes(group = hb_name), size = 1) + geom_text_repel(data = . %>% arrange(desc(month)) %>% group_by(hb_name) %>% slice(1),aes(label = round(rate_per1k, 1)),nudge_x = 0.5, size = 3, direction = "y", box.padding = 0.3,point.padding = 0.3, show.legend = FALSE) +geom_point(data = data.frame(hb_name = unique(rate_sc$hb_name)), aes(x = NA_real_, y = NA_real_, colour = hb_name), shape = 15, size = 5,show.legend = TRUE) + labs(title = "Scabies Outbreak in Scotland (June to December 2023)", subtitle = "Monthly Trends in Scabies Treatment Prescriptions", x = "Month", y = "Quatity of Antiparasitic Prescriptions /1,000 people", colour = "Health Board") +theme_bw() + guides( colour =guide_legend(title.position = "top", title.hjust = 0.5, override.aes = list(shape = 15, size = 5)))+theme(legend.position = "bottom", legend.direction = "horizontal") + scale_color_manual(values = c("#e6194b", "#ffe119", "#B0C4DE", "#f58231","#911eb4", "#46f0f0","#f032e6", "#bcf60c","#3cb44b", "#fabebe", "#0000FF", "#e6beff","#9a6324", "#606060"))
sc_squares
The analysis of preferred scabies medications aimed to identify the most commonly prescribed treatments across Scotland’s health boards from June to December 2023. The process began with isolating the relevant data fields, which included the health board, general practice, medication description, and the quantity of prescriptions for each treatment. This filtered dataset was then merged with the list of health boards to accurately label each observation with its corresponding board name. The final dataset was visualized using a bar plot that displayed the quantity of prescriptions for each type of scabies treatment across all health boards. Each health board was represented on the x-axis, while the y-axis illustrated the total quantity of prescriptions. Different antiparasitic medications were color-coded to highlight variations in treatment preferences. This visualization provided a clear and comprehensive overview of the preferred medications for scabies treatment across Scotland, revealing the distribution and trends in medication choices during the analyzed period.
separated_data <- dataset_year23 %>%imap_dfr(~ .x %>%select(hbt, gp_practice, bnf_item_description, paid_quantity) %>%filter(str_detect(bnf_item_description, scabies)) %>% mutate(month = .y))
sep_with_names_sc <- separated_data %>% full_join(nhs_board, join_by(hbt == hb)) %>%filter(!is.na(paid_quantity))
sep_with_population_sc <- sep_with_names_sc %>%full_join(population_data2023, join_by(hbt == hb))
sep_bar <- sep_with_names_sc %>% ggplot(aes(x = hb_name, y = paid_quantity))+ scale_fill_brewer(palette = "Spectral")+geom_col(aes(fill = bnf_item_description), width = .5)+theme_classic()+theme(axis.text.x = element_text(angle = 65, hjust = 1))+labs(title = "Preferred Scabies Treatments Across Scotland", subtitle = "Breakdown by health board and treatment type (June - December 2023)", x = "Health Board", y = "Number of Prescriptions", fill = "Antiparasitic Medication")
sep_bar
This bar plot illustrates that the most commonly prescribed medication for scabies was permethrin cream. However, it is important to note that over time, resistance to permethrin has been observed in scabies populations, which may affect the efficacy of the treatment, and thus more affected areas may choose alternative treatments. Additionally, the preference for different scabies treatments can vary between countries, depending on factors such as medication availability and healthcare access.
To visualize the relationship between the different scabies treatments and their distribution across months, we have created a table showing the percentage of prescriptions for each medication per health board. We will focus on the health boards that experienced the highest impact, which could depend on the rate of scabies cases in those regions.
The table below presents the proportion of prescriptions for various scabies treatments (permethrin cream, malathion lotion, ivermectin tablets, and benzyl benzoate) in June, September and December 2023. It highlights how these proportions varied across the top five most affected health boards in Scotland—Tayside, Fife, Shetland, Grampian, and Lothian— each with rates exceeding 80 antiparasitic prescriptions per 1,000 people, listed in descending order of prevalence. By calculating the percentage of each treatment prescribed relative to the total scabies prescriptions in each health board and month, we can examine trends over time. This analysis highlights the most commonly used treatments for scabies and offers insights into regional variations in medication usage.
total_presc_table <- sep_with_names_sc %>%group_by(hb_name, month) %>%summarise(total_prescriptions =sum(paid_quantity)) %>% ungroup()
presc_table <- sep_with_names_sc %>% filter(!is.na(bnf_item_description)) %>%group_by(hb_name, month, bnf_item_description)%>% summarise(treatment_prescriptions = sum(paid_quantity, na.rm = TRUE)) %>% full_join(total_presc_table, by = c("hb_name", "month"))%>% mutate(percentage =(treatment_prescriptions / total_prescriptions) * 100)%>%group_by(hb_name) %>% select(hb_name, month, bnf_item_description, percentage)%>% mutate(month = factor(month, levels = c("June", "July", "August", "September", "October", "November", "December"))) %>%filter(month %in% c("June", "September", "December"), hb_name %in% c("Shetland", "Tayside", "Lothian", "Grampian", "Fife"))%>% arrange(hb_name, month)%>%pivot_wider(names_from = bnf_item_description, values_from = percentage) %>%replace_na(list("MALATHION 0.5% AQUEOUS LIQUID" = 0, "PERMETHRIN 5% CREAM" = 0, "IVERMECTIN 3MG TABLETS" = 0, "BENZYL BENZOATE 25% APPLICATION" = 0))%>% gt()%>%tab_header(title = "Scabies Treatments by Health Board and Month", subtitle = "Data from Scottish Health Boards, June to December 2023") %>% cols_label(hb_name = "Health Board", month = "Month", "MALATHION 0.5% AQUEOUS LIQUID" = "Malathion (%)", "PERMETHRIN 5% CREAM" = "Permethrin (%)", "IVERMECTIN 3MG TABLETS" = "Ivermectin (%)", "BENZYL BENZOATE 25% APPLICATION" = "Benzyl Benzoate (%)") %>% fmt_number(columns = c("MALATHION 0.5% AQUEOUS LIQUID", "PERMETHRIN 5% CREAM", "IVERMECTIN 3MG TABLETS","BENZYL BENZOATE 25% APPLICATION"), decimals = 1) %>% cols_align(align = "center", columns = c("MALATHION 0.5% AQUEOUS LIQUID", "PERMETHRIN 5% CREAM", "IVERMECTIN 3MG TABLETS","BENZYL BENZOATE 25% APPLICATION")) %>% summary_rows(columns = c("MALATHION 0.5% AQUEOUS LIQUID", "PERMETHRIN 5% CREAM", "IVERMECTIN 3MG TABLETS","BENZYL BENZOATE 25% APPLICATION"),fns = list("Average" = ~mean(., na.rm = TRUE)),fmt = list(~ fmt_number(., decimals = 2))) %>% grand_summary_rows(columns = c("MALATHION 0.5% AQUEOUS LIQUID", "PERMETHRIN 5% CREAM", "IVERMECTIN 3MG TABLETS","BENZYL BENZOATE 25% APPLICATION"),fns = list("Overall Average" = ~mean(., na.rm = TRUE)),fmt = list(~ fmt_number(., decimals = 1))) %>%cols_move(columns = "IVERMECTIN 3MG TABLETS", after = "BENZYL BENZOATE 25% APPLICATION") %>%tab_spanner(label = "Topic Treatment", columns = c("MALATHION 0.5% AQUEOUS LIQUID", "PERMETHRIN 5% CREAM", "BENZYL BENZOATE 25% APPLICATION")) %>% tab_spanner(label = "Oral Treatment", columns = c("IVERMECTIN 3MG TABLETS")) %>% tab_style(style = cell_fill(color = alpha("#e6beff", 0.2)),locations = cells_body(rows = hb_name == "Shetland")) %>%tab_style(style = cell_fill(color = alpha("maroon", 0.2)),locations = cells_body(rows = hb_name == "Tayside"))%>% tab_style(style = cell_fill(color = alpha("lightblue", 0.2)),locations = cells_body(rows = hb_name == "Lothian")) %>%tab_style(style = cell_fill(color = alpha("#D3D3D9", 0.2)),locations = cells_body(rows = hb_name == "Grampian")) %>%tab_style(style = cell_fill(color = alpha("#B4E0B8", 0.2)),locations = cells_body(rows = hb_name == "Fife")) %>%tab_style(style = cell_fill(color = alpha("#FFF2B2", 0.7)), locations = cells_summary(columns = everything())) %>%tab_style(style = cell_text(weight = "bold"), locations = cells_summary(columns = everything())) %>%tab_style(style = cell_fill(color =alpha("#FFB84D", 0.4)),locations = cells_grand_summary(columns = everything())) %>%tab_style(style = cell_text(weight ="bold"), locations = cells_grand_summary(columns = everything()))
presc_table
| Scabies Treatments by Health Board and Month | |||||
| Data from Scottish Health Boards, June to December 2023 | |||||
| Month |
Topic Treatment
|
Oral Treatment
|
|||
|---|---|---|---|---|---|
| Malathion (%) | Permethrin (%) | Benzyl Benzoate (%) | Ivermectin (%) | ||
| Fife | |||||
| June | 41.0 | 59.0 | 0.0 | 0.0 | |
| September | 17.8 | 82.2 | 0.0 | 0.0 | |
| December | 10.3 | 89.4 | 0.0 | 0.2 | |
| Average | — | 23.05 | 76.88 | 0.00 | 0.07 |
| Grampian | |||||
| June | 47.1 | 52.8 | 0.0 | 0.1 | |
| September | 8.9 | 90.6 | 0.3 | 0.2 | |
| December | 2.6 | 97.2 | 0.0 | 0.2 | |
| Average | — | 19.53 | 80.21 | 0.11 | 0.14 |
| Lothian | |||||
| June | 41.9 | 58.1 | 0.0 | 0.0 | |
| September | 12.8 | 87.1 | 0.0 | 0.0 | |
| December | 3.1 | 96.5 | 0.1 | 0.2 | |
| Average | — | 19.27 | 80.59 | 0.04 | 0.09 |
| Shetland | |||||
| June | 69.4 | 30.6 | 0.0 | 0.0 | |
| September | 16.1 | 83.9 | 0.0 | 0.0 | |
| December | 0.0 | 100.0 | 0.0 | 0.0 | |
| Average | — | 28.52 | 71.48 | 0.00 | 0.00 |
| Tayside | |||||
| June | 25.2 | 74.7 | 0.0 | 0.0 | |
| September | 5.3 | 94.6 | 0.0 | 0.1 | |
| December | 5.6 | 94.1 | 0.0 | 0.4 | |
| Average | — | 12.06 | 87.80 | 0.00 | 0.14 |
| Overall Average | — | 20.5 | 79.4 | 0.0 | 0.1 |
This table clearly shows that permethrin cream was the preferred treatment for scabies, followed by malathion lothion. These findings are noteworthy because the World Health Organization (2023) recommends ivermectin as the primary treatment during scabies outbreaks, even in closed community settings such as hospitals, long-term care facilities, or universities. However, in the UK, ivermectin was unlicensed following the COVID-19 pandemic due to misinformation about its efficacy in treating COVID-19 and the subsequent misuse of the drug for horses (Scottish Government Health and Social Care Directorates, 2023). During the recent scabies outbreak—and likely due to increasing cases of permethrin-resistant scabies—ivermectin was relicensed in October 2023 (Scottish Government Health and Social Care Directorates, 2023). A final point to consider is that while benzyl benzoate is recognized as a potential treatment for scabies, it is unlicensed in the UK (Scottish Government Health and Social Care Directorates, 2023). This likely explains its low prescription rates.
To understand the geographic impact of the scabies outbreak in Scotland, data zones were used to localize the areas most affected within health boards (Tayside, Fife, Shetland, Grampian, and Lothian). Several challenges arose during this analysis. First, approximately 300 rows were excluded because they contained GP practice codes such as 99999, 99998, 99997, and 99996, which correspond to dentistries, hospitals, pharmacies, and unallocated codes. These entries could not be mapped to specific locations and were therefore omitted from the mapping process.Additionally, the dataset linking GP practice codes to data zones was incomplete, as not all data zones include a GP practice. The full list of data zones also lacked direct mappings to health boards, creating further challenges in localizing the impact of the outbreak. To address this issue, a supplementary dataset was used to link each data zone to its respective health board. This supplementary dataset, which was loaded earlier in the report, provided the necessary information to create accurate spatial representations. It is important to note that this analysis relied on spatial data from the 2011 data zones. While some boundary changes may have occurred since then, the supplementary dataset ensured that each data zone was accurately assigned to its current health board. These measures allowed for a detailed geographic analysis of the scabies outbreak, despite the limitations in the original datasets.
The data was prepared by merging the appropriate datasets and converting the resulting dataset into a spatial format. This spatial data was then transformed to a latitude and longitude coordinate reference system (CRS) to ensure compatibility with the leaflet mapping package. A tibble was created to define health board codes, names, and their geographic coordinates. To enhance visualization, a function was implemented to create a color palette based on the logarithmic transformation of the paid_quantity variable (log1p()), reducing skewness and highlighting variations across data zones. Another function was developed to generate an interactive map for each specified health board. Each map was centered on the health board’s geographic coordinates, and the legend was back-transformed to the original scale using expm1() for clarity. Finally, the maps were combined and formatted into a grid layout.
with_datangp <- total_data23_sc %>% full_join(gp_scot, join_by(gp_practice == practice_code))%>% full_join(data_zones, join_by(data_zone == DataZone))%>% full_join(DataZone_with_hb, join_by(data_zone == dz2011_code)) %>% mutate(paid_quantity = replace_na(paid_quantity, 0)) %>% select(-hbt, -hb)
with_datangp_sf <- st_as_sf(with_datangp)
hb_data <- tibble(hb_code = c("S08000026", "S08000030", "S08000024", "S08000020", "S08000029"),title = c("NHS Shetland", "NHS Tayside", "NHS Lothian", "NHS Grampian", "NHS Fife"),lat = c(60.7931, 56.4852, 55.9444, 57.1497, 56.0656),lng = c(-1.3309, -3.0653, -3.1785, -2.1095, -3.1472))
with_datangp_sf <- st_transform(with_datangp_sf, crs = 4326)
create_color_palette <- function(data) { log_paid_quantity <- log1p(data$paid_quantity)
colorNumeric(palette = "Blues",domain = log_paid_quantity,na.color = "transparent")}
create_leaflet_map <- function(hb_code, title, data, lat, lng) {
filtered_data <- data %>% filter(hb_code == !!hb_code)
color_palette <- create_color_palette(filtered_data)
leaflet(filtered_data) %>%addProviderTiles("CartoDB.Positron") %>% addPolygons(fillColor = ~color_palette(log1p(paid_quantity)),color = "white", weight = 1, opacity = 1, fillOpacity = 0.8,popup = ~paste("<strong>", title, "</strong><br>", "Paid Quantity: ", paid_quantity)) %>% addLegend(position = "bottomright",pal = color_palette, values = log1p(filtered_data$paid_quantity),title = "Quantity",labFormat = labelFormat(transform = expm1, digits = 0),opacity = 0.5,labels = ~round(expm1(log1p(filtered_data$paid_quantity)), 0)) %>%setView(lng = lng, lat = lat, zoom = 8)}
maps <- purrr::map2(hb_data$hb_code,hb_data$title, ~create_leaflet_map(.x, .y, with_datangp_sf, hb_data$lat[hb_data$hb_code == .x], hb_data$lng[hb_data$hb_code == .x]))
combined_map <- div(style = "display: flex; flex-direction: row; justify-content: space-between; flex-wrap: wrap; gap: 10px;", lapply(1:length(maps), function(i) { div(style = "margin-right: 20px; width: 250px;",h4(hb_data$title[i]),maps[[i]])}))
browsable(combined_map)